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A finite quantum system evolving unitarily equilibrates in a probabilistic fashion. In the general many-body 
setting the time-fluctuations of an observable A are typically exponentially small in the system size. We consider 
here quasi-free Fermi systems where the Hamiltonian and observables are quadratic in the Fermi operators. We 
first prove a novel bound on the temporal fluctuations and then map the equilibration dynamics to a generalized 
classical XY model in the infinite temperature limit. Using this insight we conjecture that, in most cases, a 
central limit theorem can be formulated leading to what we call Gaussian equilibration: observables display a 
Gaussian distribution with relative error AA/A = 0(L -1 / 2 ) where L is the dimension of the single particle 
space. We prove the conjecture analytically for the magnetization in the quantum XY model, and numerically 
for a class of observables in a tight-binding model. 
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Out of equilibrium and equilibration dynamics of closed 
quantum systems have been recently at the center of a renewed 
and intense interest OUT). In particular the issue whether 
quantum integrability plays a key role in equilibration and, if 
so, which one, has been investigated by several authors ||8}- 
Kfl . In this paper we will address this problem along the 
lines of the probabilistic approach to quantum equilibration 
advocated in 171 [Til . Here the central object is the (infinite 
time) full time statistics of the expectation value of a quan- 
tum observable. We will focus on a particular, yet very impor- 
tant class of quantum integrable systems: quasi-free fermionic 
systems. Using fairly general central-limit type arguments as 
well as explicit analytic examples we will argue that a sharp 
distinctive feature of these systems, as opposed to the gen- 
eral interacting ones, is an exponential enhancement of the 
amplitude of the temporal fluctuations of a quadratic observ- 
able around its mean value. This appears to be a precise and 
quantitative way to make sense of the common folklore that 
integrability leads to a poorer (or no) equilibration. 

The system is initialized in the generic state p with N 
particles and both the evolution Hamiltonian H and the ob- 
servable A are quadratic in the fermionic operators EUll . The 
Hamiltonian is H = Y^xy c xM x ,yC y — c^Mc (notation 

c< = ^c{, . . . , cQ, L number of sites). The general quadratic 

observable has the form A = J2 X v c t 

We 

will assume that WaW^ = O (1) ED as this guarantees 
that the expectation values of A scales at most extensively 
with the system size [22 1 The main object of investigation 



is A (t) = tr {Ae~ ltH poe ltH ) . Exploiting the Gaussian na- 
ture of the problem and introducing the covariance matrix 
R y x := tr {poc\.c y ) (0 < R < J) one can show that the 
expectation value A (t) reduces to a trace in the one-particle 
space: 



A (t) = tr (ae 



-itM 



Re 



(1) 



Eq. ([T} is perfectly analogous to its many-body version 
A (t) = tr {Ae~ itH p e ltH ) with R playing the role of the 



initial state po- There is however one importance difference: 
while trpo = 1 one has trR = N = vL , i.e. is extensive (we 
defined v = N/L the filling factor). 

For unitary evolution, the density matrix p (t) = 
e~ ltH p a e ltH does not converge neither in the strong nor in the 
weak topology |6|. Equilibration must be formulated in prob- 
abilistic terms. Given the observation time window [0, T], the 
observable A (t) has probability Pa {a) da of being in the in- 
terval [a, a + da], where the probability density is given by 

Pa {at) — S {A (i) — a) T , and the time average operation is 
jt — rp-i jT j (j-^ ^ p or s i m pii c ity we will always take the 
limit T — » oo when taking time average and write simply / in 
place of f°°. Roughly speaking an observable A equilibrates 
if its probability density Pa {a) is highly peaked around its 
mean A {t) — tr {A~p) . 0. The role of equilibrium state is 
played by the time-averaged density matrix ~p. 

The question we are going to address here is: what is the 
size of the fluctuations of Pa {a) for observables A in this 
quasi-free setting? 

Let us start with a few remarks on the nature of the equilib- 
rium state itself. If the spectrum is non-degenerate, the aver- 
age, dephased, state p has the form ~p = J2nPn\ n )( n \ where 
\n) are many -body eigenstate of H corresponding to energy 
E n and p n = (n\po\n). The powers of the Hamiltonian H n , 
n = 0, 1, . . . , d— 1 are linearly independent if there are d dif- 
ferent eigenvalues E n and if the spectrum is non-degenerate 
d coincides with the Hilbert space dimension [23 1. In this lat- 
ter case the average state can always be written in the form 

p = cxp X)fe=o a kH k , the so called generalized Gibbs en- 
semble [81. The coefficients a n depend on the initial state po 
and on the eigenvectors {|n)}. The condition to write the co- 
efficients a n in terms of the p n is precisely that of the invert- 
ibility of the Vandermonde matrix V n ,k = {En)'' -1 (n, k = 
1, ... ,d), i.e., once again, non-degeneracy of the spectrum, 
since det V — Iln<fe i^k ~ E n ). The relation expressing the 



p n 's in terms of the a^'s is: p n = exp 



Since 
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P = expX;„lnp,»(n| = expY,k=aT, n a k E %\ n )(n\ ' the 
inverse relation is, in vector notation a = V^ 1 In p. The in- 
verse of the Vandermonde matrix can be found explicitly by 
expanding the identity |n)(n| = Ilfe^n E~-^E k ■> multiplying 
by \np n , summing over n and exponentiating see [24 1 

Let us now go back to the context of Eq. ([TJi. The time av- 
eraged covariance R — e~ ttM Re ltM defines a Gaussian state 
pg, in fact < R (t) < 1 implies < R < I. Moreover, 
since tri? (<) = N, for all t, tri? = N, p-^ is a Gaussian state 
with N particles. Now, for what concerns quadratic observ- 
ables of the kind A = c) ac, their time average expectation 
value is the same as that obtained with p^\ A — traR — 
tr Ap-ft . In other words the states p and p-^ are the same when 
restricted to quadratic observables. A generic Gaussian state 
can be written in the form p R — A/"exp ^0^)6) where 
G{R) := log[i? T (l - R T )~ 1 } [12] and TV is a normalization 
constant. From this it immediately follows that that p-^ can 
be written as p-g = /VexpY^ \ k c\c k where c k s are eigen- 
modes of H and coincides with ppc of (S). When the system 
equilibrates i.e., A(t) — > tr (Ap^) this remark shows the va- 
lidity of the GGE (in the sense of (8)) for any initial state po 
and all quadratic observables. 

A bound on the variance Assuming the non-resonant con- 
dition (E n — E m = E p — E„ implies n — m and p = q 
or n — p and m — q), Reimann has shown iTPTl that the 

temporal fluctuations AA 2 = (A (i) - A) 2 satisfy AA 2 < 
diam (A) 2 tip 2 where diam (A) is the maximum minus the 
minimum eigenvalue of A. Now typically, for most initial 
states po, the purity trp 2 is exponentially small in the system 
size implying exponentially small fluctuations 0. 

In the quasi-free setting the non-resonant condition does 
not hold. Let us then seek for the analogous of the 
bound of Reimann in our quasi-free case. Let the one- 
particle Hamiltonian have the following diagonal form M — 
J2k A-k\k)(k\. The time averaged covariance matrix is 

then R = T,k( k \ R \ k )\ k )( k l We a l so deflne F k,g = 
(k\a\q) (q\R\k) . Assuming the non-resonance of the gaps 
of the one-particle spectrum, one gets AA 2 = trF 2 — 
E k (Fk, k f < trF 2 = Z k J(k\a\q)\ 2 \(q\R\k}\ 2 . Now 
R is a non-negative operator and so induces a (possibly 
degenerate) scalar product which satisfies Cauchy-Schwarz 
inequality: |(<?|i?|fc)| 2 = \(q\k) R \ 2 < (q\q) R (k\k) R = 
(q\R\q)(k\R\k). This leads us to 

AA 2 < tr (aRaR) < Waf^ tri? 2 (2) 

Now, since < R < 1, tr^R 2 < tr^R = tri? = N, we 
finally obtain AA 2 < \\a\\ 2 xi i'L . While equation (j^Ji is the 
quasi-free analog of the Reimann's bound it implies some 
important differences with respect to the general (non-free) 
case. Consider the situation where the observable A is ex- 
tensive. In the non-free case the diameter of A is extensive, 
i.e. diam (A) = O (L D ^j in D spatial dimensions. More- 
over, the minimum value of the purity tr~p 2 is 1/d and so 
is exponentially small in the system size. In the quasi-free 



setting, instead, the minimum value of tri? in Eq. &2\ is 

2 — 

mintri? = N 2 / L = v 2 L, since R/N defines a density ma- 
trix for which the minimum purity is the inverse of the matrix 
dimension L. All in all, recalling that WaW^ = O (1) , for 
extensive observables, the minimum of the bound to the vari- 
ance is O (L 2 e~ 7L ), 7 > 0, (for systems of linear size L), 
whereas in the quasi-free setting one has O (L). This seems 
to hint at the fact that fluctuations in the quasi-free setting are 
proportional to the system size and are hence much larger than 
in the non-free case where they are exponentially small in the 
volume. Of course Eq. |2]) is just an upper bound, and nothing 
prevents, in principle, from having a much smaller variance. 
In the following, however, we will argue that the extensive 
behavior of the fluctuations, in quasi-free systems, is in fact 
quite general. 

Mapping to a classical XY model Let us write again the 
generic expectation value ([T]) in the basis which diagonalizes 
M: A (t) = A+2 J2 k<q \F k , q \ cos (t (A k - A q ) + (f, k>q ) with 
<t> k .q — ar g-Ffe,9- To obtain information on the probability 
density Pa (a) we consider the generating function \A (A) :— 
e A (-4(*)--4) Now we observe that if the (one-particle) ener- 
gies are rationally independent (RI), as a consequence of a 
theorem on the averages, the infinite time average of A (t) is 
the same as the uniform average over the torus T L . In this 
case the generating function \A (A) is exactly given by the 
partition function of the generalized, classical XY model with 
energy E (i?) = 2 Y. k<q \ F k,q\ cos (i9 fc - d q + M ) and in- 
verse temperature f3 = —A. The matrix \F k!q \ defines the 
lattice of the interactions while the phases <f> k:q give the offset 
from which the angles are measured. Note that the behavior 
of the density Pa (a) is dictated by \A (A) in a neighborhood 
of A = which corresponds to infinite temperature of the 
classical XY model. 

It is not difficult to engineer a situation which exactly repro- 
duces the standard XY model in £>-dimension. For example, it 
suffices to consider the Hamiltonian H = J2 X li>xcl.c x , with 
p x RI (x is a point of a £>-dimensional lattice), choose the 
observable A = ^ c\,Cy ((x, y) indicates nearest neigh- 
bor) and initial (Gaussian) state IV'o) = i -1 ^ 2 l~2 x \x) = 
c jt=ol0)- In this case the partition function Z (free-energy 
J 7 ) of the classical £>-dimensional XY model is precisely the 
characteristic function of the observable A: \a (A) = Z 
(\?tXA (A) = T). In fact in this case a x y = 5^ x y ^ while 
Rx.y = 1/L so that F x y = L~ 1 8( x ^ y ) which defines the 
nearest neighbor hyper-cubic graph. 

We would like to stress here that the one-particle space has 
a natural underlying geometric structure. For instance, the 
labels k, q represent points in momentum (real) space in a 
superfluid (localized) phase and the distance |fc — q\ is well 
defined. Now, when the matrix elements |i*fc g | decay suffi- 
ciently fast as \k — q\ — > oo the corresponding XY model is 
well defined in the thermodynamic limit, i.e. the intensive free 
energy has a limit as L — >• oo. This happens for instance in 
case |-Ffc.q| decays decays exponentially in |fc — q\ or if one 
has l-Ffc^l ~ 1/ |fc — g| 7 with 7 > D. When this is the case 
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one has \A (A) = exp L D T (A). Moreover one expects T (A) 
to be analytic in the high temperature, A = 0, limit, meaning 
that all the cumulants are extensive. From this we immedi- 
ately draw the central limit theorem (CLT): as L — > oo the 
variable (A (t) — A) / L D / 2 tends in distribution to a Gaussian 
with zero mean and finite variance given by d\_ F (X). We 
call this situation Gaussian equilibration. It is important to 
stress that one cannot have Gaussian equilibration in the non- 
free setting or all the cumulants would scale in the same way. 
Instead in the non-free setting one expects extensive average 
but variance exponentially small in the volume. 

We will now further corroborate these arguments with two 
examples where indeed Gaussian equilibration can be proven 
or shown. 

Quench on the quantum XY model The Hamiltonian is 
given in terms of Pauli spin operators a x ' v ' x (we use periodic 
boundary conditions) 



H 



L 

£ 

i=l 



1 + 7 



1 



'7 



'i+l 



hat 



(3) 

In the quench scenario, the initial state is the ground state 
\i{jo) of the Hamiltonian with certain values of the parameters 
7o, h . The parameters are then suddenly changed and |^ ) is 
evolved with the Hamiltonian with parameters 71, hi. 

The model in Eq. <j3j has been long used as a testbed for the 
study of quantum phase transitions (taking place at h = 1 and 
zero temperature) in many body systems, and more recently in 
the realm of out-of equilibrium unitary dynamics. See lfl4l[T5l 
and the more recent monumental [ 16 1 for more details and ref- 
erences. A Jordan- Wigner transformation brings Eq. ^ to a 
quadratic form in Fermi operators. Since of in terms of Fermi 
operators is of = 2c\ci — 1, the transverse, total magnetiza- 
tion M(t) = Y^ii^i (t)) i s a quadratic observable. Its expec- 
tation value in the quench setting is given by lfT31 : M (t) — 

2 E fc >o cos ^fe 



(i) 

k 



sin'tfjj. sin (S-dk) 



cos (S'&k) 

= — 7^ sinfc/ (hi + cos k), d'&k 
\J (7< sinfc) 2 + (hi 



K) 

«<■> - 



where tan$ n 

•d^ and = 2\J (7^ sinfc) z + (hi + cos/c) z are the one- 
particle energies. The quasi-momenta are quantized according 
to k = 7T (2n + 1) /L, n = 0, 1, . . . ,L/2 - 1. At this point it 
seems quite natural to expect that the energies A^ are ratio- 
nally independent. For example one can show that the num- 
bers cos(k n ) k n — tt (2n + 1) /L, n = 1, 2, . . . , (L — l)/2 
are rationally independent for L prime Q. Given the form of 
the dispersion we may expect that the requirement that L is 
prime may be lifted. 

Assuming rational independence of the one-particle ener- 
gies A[^, the corresponding classical XY model has energy 

E W = Efc>o^ cos (^fe) with = sin4 1) sm(W fe ). 
Each classical spin ' k' interacts with an external field along a 
fixed axis with strength W k - The partition function factorizes, 
each integral over dk giving a Bessel function Iq (XWh) and 

exp E fc >o ln ( J o (AW fc )). Clearly 



M (t) — M is a sum of independent random variables, each 
with zero mean and variance /2. Now we have the follow- 
ing 

Theorem For any value of parameters (X ,h ) 7^ 
(Xi,hi), the variable (M — M)/\L as L — > 00 tends 
in distribution to a Gaussian with zero mean and variance 
(27r)- 1 J^W^/2)dk. 

We can prove the theorem by showing that the Lya- 
punov condition is satisfied so that the central limit theo- 
rem follows from Lindeberg's theorem (see e.g. ifTTl ). Fol- 
lowing the notation of [17] we have s\ — J2k a k = 
E k>0 W j/2 L(4ny 1 S^Wldk. Then, with 6 = 1 

|Afc| 2+(5 = 4/(37r) I Vt^fc 1 3 . The Lyapunov's condition with 
5 = 1 is the vanishing of the following quantity as L — >• 00: 

w (4) 
^^0 (l^\W k \ 2 dk) 



1 16 fo\W k \ 3 dk 



we obtain e x ( M ^~ M ) 



Indeed the RHS of Eq. |4]) goes to zero as L — >• 00 since 
\W k \ < 2 for all k and \W k \ ^ for almost any k for 
(Xo,h )^(Xi,hi). 

Remark For small quench close to a critical point, the 
function W k — sin('d^)(d'dk/dx) dx (x is the quenched 
variable, x = X, h), becomes highly peaked (for instance, the 
peak is around k — tt close to the Ising critical point at h = 1 
and Wk/dx diverges as 1/fc). For finite L and sufficiently 
small quench, few terms W k dominate and one can obtain (for 
finite L) a non-Gaussian distribution. Indeed, as discussed in 
detail in IfTTl . this is the case in general: for small quenches 
close to a quantum critical point, observables become a sum 
of few independent random variables and the distribution ac- 
quires a universal double-peaked form. 

Tight binding model We consider here the 1 -D tight bind- 
ing model H = J2 x (c) x c x+ i + h.c.) with twisted bound- 
ary conditions cl+i — c\e %BL as proposed in lfT8ll . As 
quadratic observable take A with diagonal one-particle ma- 
trix: A = Y^,x=i c ! c x- The system is initialized setting all the 
N particles say to the left of the chain, i.e. the initial covari- 
ance matrix is R — diag (1, 1, . . . , 1, 0, . . . , 0) with ones 
on the diagonal. The observable A is extensive for I = aL, 
and the thermodynamic limit is given by a = l/L, v = N/L 
constant and L — > 00. The time evolved observable reads 
Mt) = Ek, q 9N(k-q)ge(q-k)e- it (^-^ where the 

function g e (£) is given by g e (f) = L^ 1 Y? x =i e ~ rxi - The 
matrix F kiQ depends only on the difference k ~ q: Fk. q = 
f (k — q) = (k — q) ge. (q — k). The eigen-energies are 
given by Afe = 2 cos (k + 0) and the quasimomenta can be 
considered quantized as k = 2im/L, n = 0,l,...,L — 1. 
For 8 = the energies are degenerate as A^ = A^-k but for 
most of the 9, tt — k does not belong to the Brillouin zone and 
the energies are non-degenerate. In this case the average is 
A = Ylk 9 N (0) 9 l (0) = Lva. Still the A^ cannot be RI as 
= trM = J^k ^-k> however it is likely that there are few re- 
lations among the energies. For instance Tasaki showed lfT8l 
that for L odd and most of the 6 e (0, tc/L), the one-particle 
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energies satisfy the non-resonant condition. This implies that, 
for most 6, the variance is given by 

A.4 2 = \9N (k ~ q)gt (q ~ k)\ 2 - £ [g N (0) g t (0)] 2 . 

k,q k 

The double sum above can be evaluated going back to real 
space with result L~ 2 [Ni 2 - £ {l 2 - l) /3] (assuming N > 

1 otherwise swap N with £), hence also the variance is exten- 
sive and one has AA 2 = L [va 2 - ^ - v 2 a 2 ^j + O 

On the other hand, using methods borrowed from statisti- 
cal mechanics, we can get an approximate analytic expres- 
sion for the whole cumulant generating function. We first 
assume rational independence of the energies Afc so that the 
problem is mapped to an classical XY model. A single 
(or a few) relation among the energies is a sort of bound- 
ary condition for the classical model and is not expected 
to change the leading, bulk, term. The energy of the clas- 
sical model for the translation invariant case is E (i?) = 

2EitT 1)/2 |/(^)IELiCOs((^ +d -^) {<j> d = 

arg/ (d) and <& x is a periodic extension of d x , i.e. d x +nL 
• ' Now we note that / (d) is highly peaked 
around d = 0, so we approximate the energy keep- 
ing only the nearest neighbor term, i.e. E (i?) f=s 

2 1/ (1)1 J2x=i cos (J&x+i - &x) + 4>i)- Tms is precisely 
a one-dimensional (classical) XY model with periodic bound- 
ary condition (and off-set <f>{) and can be solved via trans- 
fer matrix method |[T9l . The partition function becomes Z = 
tr/C L where the transfer matrix operator is [K.h('& x )]('& x ^.i) = 
J o 2 "^e A|/(1)|cos (^+ 1 -^ +01 )/ l (^)/(2 7 r). The transfer 
matrix is non-hermitian because of <f>\. Using the identity 
e Kcos(a-/3) _ J2 peZ I p (K) e lp( ~ a ~^ where I p are Bessel 
functions (they satisfy Iq > 1% > I2 ■■ ■ and I p (K) — 
I- p (K)) one sees that plane waves e im,&x are eigenfunc- 
tion of £ with eigenvalue I m (2 \f (1)| A) e - im<i> . The largest 
eigenvalue in modulus, with m — 0, gives immediately the 
free energytJ 7 = hxZ/L = In (Iq(X |/(1)|). As expected, 
in this approximation, the cumulant generating function is ex- 
tensive and analytic in A = 0. One then has again CLT in its 
standard form: the variable [A (t) — A)/ \/L, tends in distri- 
bution to a Gaussian in the thermodynamic limit with variance 
d 2 =0 J- = |/(1) | 2 /2. The Gaussian prediction is clearly con- 
firmed by a numerical experiment see fig.[T] 

Conclusions In this paper we addressed the question of 
equilibration in quasi-free Fermi systems. The initial state is 
a general state while the evolution Hamiltonian as well as the 
observable are quadratic in the Fermi operators. While for 
general non-free systems the variance AA 2 is typically expo- 
nentially small in the volume, we find that in the quasi-free 
setting (for extensive observables) both the mean and the vari- 
ance are proportional to the volume. This hints at the possibil- 
ity of an underlying central limit theorem, a circumstance that 
we termed Gaussian equilibration. In this case the properly 
rescaled observable becomes Gaussian in the large L limit, 
and in general the relative error satisfies VAA 2 /A ~ 1/VT. 
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Figure 1: Full distribution for the observable A = X^=i c x c x 
around its mean for L — 201, N = 100, I = 101. The histogram 
is obtained sampling A (t) at 240,000 random times uniformly dis- 
tributed in [0, T max ] with T max = 180, 000. The thick curve is a 
Gaussian with zero mean and variance a 2 — 0.0208L as computed 
in the text. 

We prove Gaussian equilibration for the magnetization in the 
quantum XY model and give evidence for a class of observ- 
ables/initial states evolving with a tight-binding model. The 
key insight is a mapping of the equilibration dynamics to a 
generalized classical XY model at infinite temperature. As a 
by product we obtain a quantum setting (initial state, Hamil- 
tonian, and observable Axy), such that the equilibration dy- 
namics of Axy gives the solution of a classical XY model in 
D-dimension, and vice-versa. 

Ultimately the origin of the exponential decrease of the sig- 
nal to noise ratio illustrated in this paper lies in the exponential 
reduction of the effective phase space entailed by the quasi- 
free nature of our setup i.e., the many body space gets effec- 
tively replaced by the single particle one. A natural question 
for future research is whether a similar mechanism is at place 
for a more general class of integrable systems. For example 
it is tempting to conjecture that this might be indeed the case 
for systems integrable by Bethe Ansatz as interacting eigen- 
states are described in terms of single-particle quasi momenta 
(rapidities). 

Finally, we would like to stress that the enhancement of the 
temporal fluctuations of AA 2 (for a quadratic observable A) 
should provide a universal and experimentally testable signa- 
ture of integrability in the context of out of equilibrium dy- 
namics. 
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